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The existence and tlie properties of self-bound quark matter in the NJL model at zero tem- 
perature are investigated in mean-field approximation, focusing on inhomogeneous structures 
with one-dimensional spatial modulations. It is found that the most stable homogeneous so- 
lutions which have previously been interpreted as schematic quark droplets are unstable 
against formation of a one-dimensional soliton-antisoliton lattice. The solitons repel each 
other, so that the minimal energy per quark is realized in the single-soliton limit. The prop- 
erties of the solitons and their interactions are discussed in detail, and the effect of vector 
interactions is estimated. The results may be relevant for the dynamics of expanding quark 
matter. 
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I. INTRODUCTION 

The Nambu-Jona-Lasinio (NJL) model [])] is a popular tool for studying low-energy properties 
of strongly interacting matter, like spectra and scattering of light hadrons, or the phase diagram 
at nonvanishing temperatures or densities (for reviews, see Refs. 0-0] )• While being relatively 
simple, the NJL model shares the global symmetries of QCD, in particular chiral symmetry, which 
is considered to be the most important feature of the model. On the other hand, it is well known 
that the NJL model lacks confinement. In this sense it can be viewed as complementary to the 
MIT bag model 0], which is confining by construction, but violates chiral symmetry at the surface. 

Some time ago, it was realized, however, that for sufficiently attractive interactions, the NJL 
model at zero temperature has solutions of self-bound chirally restored quark matter, which can 
be interpreted as bag-model-like quark droplets 0, 0-0] • In fact, the link between both models 
is the existence of a "bag pressure", which in the bag model is introduced by hand in order to 
stabilize the solutions, whereas in the NJL model it is a dynamical consequence of spontaneous 
chiral symmetry breaking in vacuum. 

The self-bound quark matter solutions mentioned above have been obtained in the thermody- 
namic limit and correspond to infinite homogeneous matter. Their interpretation as quark-matter 
droplets is based on the behavior of the energy per particle, E/N, which shows a minimum at 
some nonvanishing saturation density. This means, a finite piece of quark matter with this density 
would be stable against collapse or expansion, just like a liquid drop. 

An equivalent statement is that the matter has vanishing pressure, which means, it is in me- 
chanical equilibrium with the vacuum. Thus, in order to have a solution of this type, there must 
be a phase coexistence of the vacuum with a dense-matter phase. In other words, at some critical 
chemical potential, there must be a first-order phase transition from the vacuum to dense matter. 
This is realized in the NJL model, if the interaction is sufficiently attractive. On the other hand, 
if the attraction is relatively weak (a condition which can be achieved, e.g., by adding a repulsive 
vector interaction), it is also possible to have a second-order phase transition or a crossover. In this 
case there is no stable matter solution and E/N takes its minimal value at zero density. Without 
applying external forces, a finite piece of quark matter would then keep expanding, i.e., behave like 
a gas. 

It is tempting at this point to extrapolate the self-bound solutions down to droplets consisting 
of only three quarks, and to interpret them as baryons. However, although some of the resulting 
"baryon" properties are quite reasonable [7|, it is obvious that this extrapolation is not reliable. In 
fact, using solutions for infinite homogeneous quark matter to describe finite droplets, one has to 
assume that surface effects can be neglected. This assumption might be justified for large droplets 
but most likely not for small ones. Besides, if the surface tension is positive, as derived, e.g., in 
Refs. 1CH12||. smaller droplets are disfavored. The preferred state in the model would therefore 
be a configuration where all quarks are joined in one big spherical nugget, rather than hadronized 
into individual baryons. 

It turns out, however, that this is not quite the case. More recent studies of the NJL phase 
diagram have revealed that the first-order chiral phase transition between homogeneous phases 
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gets replace d by an inhomogeneous region if one allows the chiral condensate to be nonuniform 
in space 13l4l5||. In this region, a special class of solutions, which vary in one spatial dimension, 



has been found to be favored over all other shapes considered so far. These solutions correspond 
to a lattice of domain-wall soliton^ described in terms of Jacobi elliptic functions and smoothly 
interpolate between the homogeneous chirally broken and restored phases [l5| . In particular at 
the low-density side, they take the form of a single soliton, which is thermodynamically degenerate 
with the homogeneous chirally broken phase. As a consequence, the phase transition to the latter 
is of second order. 

This changes our picture of self-bound quark matter in the NJL model considerably. Since 
there is no longer a first-order phase transition connecting the vacuum with a finite-density phase, 
but a second-order phase transition to the inhomogeneous phase, the minimal E/N should now 
be reached at zero average density. However, unlike the homogeneous case where the low-density 
regime corresponds to a dilute gas of constituent quarks, we now expect a "liquid crystal" of well 
separated solitons and anti-solitons. These objects have a nonvanishing quark density and a finite 
size in one spatial dimension, while being infinite in the remaining two dimensions. This could be 
seen as a step towards "real" quark droplets, which are finite in three dimensions. 

In the present article we perform an explicit model study to investigate the scenario outlined 
above quantitatively. After briefly introducing the formal background, we calculate E/N as a 
function of the average density and compare the results for inhomogeneous solutions with those 
for homogeneous matter. Based on these results, we then discuss the properties of single solitons 
and their interactions. Finally, we estimate the effect of vector interactions, before we draw our 
conclusions. 



II. NONUNIFORM QUARK MATTER IN THE NJL MODEL 

In this section we briefly summarize the main properties of the one-dimensional solitonic NJL- 
model solutions derived in Refs. and 0. Afterwards we study the single-soliton limit of these 
expressions. 

A. Mass functions and thermodynamic potential 

Our starting point is the Nambu-Jona Lasinio Lagrangian [1] in the chiral limit, 

Cnjl = i^il^d^i^ + G [{^-^f + {^i^^'Ta^f) , (1) 

where -0 denotes a quark field with Nf = 2 flavor and Nc = 'i color degrees of freedom, Ta are the 
Pauli matrices in flavor space, and G is a dimensionful coupling constant. The model is studied 
in the mean-field approximation. To this end, we assume the presence of a nonvanishing scalar 
condensate, {ipip) = S{z), which we allow to vary in one spatial dimension (z direction) while 



^ For brevity, we will just call them "solitons" in the following. 
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being constant in the two perpendicular directions {x and y) and in timel§ Accordingly, the quarks 
acquire a z-dependent dynamical mass function M(z) = —2GS{z). 

With this ansatz, one can employ the known results for the 1+1-dimensional Gross-Neveu 
model 1^, to construct solutions of the 3+1-dimensional problem [l^. The mass function can be 
expressed in terms of Jacobi elliptic functions, 



M{z) 



sn(Az|i^)cn(A2;|i/) 



(2) 



dn(Az|z^) 

characterized by two parameters: an amplitude A and the so-called elliptic modulus G [0, 1]. The 
latter determines the shape of the modulation, continuously changing from sinusoidal for = to 
a hyperbolic tangent ("kink") for v = 1. For u < 1, M{z) is periodic with period [Ji 
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where K is the complete elliptic integral of 1st kind. 

For the thermodynamic potential per volume at temperature T and chemical potential ^ one 
obtains 



Q(T, ^l■ A, u) = -NfN, / dE p,nh{E; A, u) [Uac{E) + fmed{E- T, ^)] 
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Here K is again the complete elliptic integral of 1st kind, F is the incomplete elliptic integral of 
1st kind, and E are the (complete or incomplete) elliptic integrals of 2nd kind. Furthermore we 
introduced the notations 9 = 1 — u, 9 = arcsin(ii^/(-v/pA)), and 9 = arcsin(A/ii^). 
The functions fvac and fmed{E) in Eq. (jl]) are given by 



and 



fvac{E) — E 



/™ed(^;T,//) =Tlog( l + exp( ) ) + riog( 1 + exp' 



(6) 



(7) 



Other cases, like chiral density waves, which also include pseudoscalar condensates [l^. [l^. or two-dimensional 
crystals [l^l have been considered as well, but have been found to be less favored at low densities [isl. [itI. [isll . 



5 



> 

0) 




2.5 
2 

1.5 
1 

0.5 




10 15 20 
z(fm) 



30 



10 15 20 25 30 
z(fm) 



FIG. 1. Left: Mass function M{z) for A = 330 MeV and ^ I - IQ-^^. Right: Corresponding density 
profile at T = and ^ = 323.3 MeV. 

Since the vacuum part of the energv integral is divergent, we have to regularize it. We use Pauli- 
Villars regularization of the form [3|] 



j=0 



(8) 



with Co = 1, ci = —3, C2 = 3, C3 = —1 and a cutoff parameter A. 

With these expressions at hand, the ground state of the system can be determined by minimizing 
the thermodynamic potential in the two parameters A and z^. 



B. Density profile 

The density profiles of the above solutions are given by [1 



niz) = NfN,j dEpD,inhiE,z;A,i^){n+{E)-n^{E)) , (9) 



where 

are the Fermi occupation functions for particles and antiparticles, respectively, and the density 
matrix element PD,inh can be related to Pinh, Eq. ([5]), upon the replacement 



As an example, we show in Fig. [T] the mass function for A = 330 MeV and z/ = 1 — 10^^^ (left), 
and the corresponding density profile at T = and p = 323.3 MeV (right). Comparing these 
figures, one can see that the density is peaked at the points where the mass functions vanish, i.e.. 
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the regions of high density correspond to the regions where chiral symmetry is almost restored. 
This is reminiscent of the bag model, where the quarks are only allowed in the trivial vacuum. 

Because of the alternating sign of the mass function, the density peaks can be identified with 
solitons and antisolitons. As obvious from Eq. pT]) . the density does not depend on the sign of the 
mass function. Therefore the distance a between two neighboring peaks is equal to one half of the 
period L, where L is given in Eq. ([3]), 

(-) 

C. Single-soliton limit 

In the limit i/ — )• 1, the period L goes to infinity, and the mass function, Eq. features a 
single kink at z = 0, 

M{z)\u=i = A tanh(Az) , (13) 

corresponding to a single soliton. The density of states, Eq. ([5]), becomes 

p,nh{E; A, 1/ = 1) = e{E - A) -^eVe^ - = phom{E; A) , (14) 

which is equal to the density of states of an ideal gas of quarks with constant mass A. As a conse- 
quence, the free energy of the inhomogeneous phase in the single-soliton limit becomes degenerate 
with the free energy of homogeneous matter with a constituent quark mass A. 
For the density matrix element Eq. (|lip one gets 



PD,mh{E, z;A,u = 1) = phom{E; A) + pD,soi{E, z; A) , (15) 

with a localized part 

PD,soi{E, z;A) = ^ (e{A -E) + 9{E - A) - arcsin ^ ) ] (16) 
47r \ TT E J cosh (Az) 

and a homogeneous background given by Eq. (jl4[) . The latter is again equal to the analogous term 
in homogeneous matter. 

Accordingly, the density one obtains from Eq. ([9|) can be separated into a constant background, 
which is equal to the density in a homogeneous ideal gas of quarks with mass A at given temperature 
and chemical potential, and a localized peak, which corresponds to the extra quarks in the solitons. 
In particular, since the localized part drops off exponentially at large values of \z\, the average 
density 

L 

n = lim — / dz n{z) (17) 



-L 



is entirely determined by the background and, thus, equal to the density in homogeneous matter. 
As a consequence, a phase transition from the inhomogeneous phase to the homogeneous chirally 
broken phase taking place at = 1 is second order. 
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The additional density contribution due to the quarks in the sohtons, 

oo 

nsoi{z) = NfN, J dEpD,soi{E,z-/^){n+{E)-n,{E)) , (18) 



is perhaps the most interesting part. In particular, nsoi{z) and, thus, n{z) is nonzero even at T = 
and ^ < A, when the background density vanishes. In this case, which corresponds to a single 
soliton embedded in the vacuum, one finds 

n.oK^)kA>.,T=o} = ^^^^^^. (19) 

Note, however, that here we have simply assumed that solutions with v = 1 and A > /i exist at 
zero temperature. Of course, we have to check whether this comes out of the minimization of the 
thermodynamic potential. Since = 1 is realized exactly at the second-order phase transition 
from the homogeneous to the inhomogeneous chirally broken phase, this means that at the critical 
chemical potential Hc,inh, the amplitude A must be bigger than fic,inh- On the other hand, at 
A* = tJ-c,inh, the amplitude A is equal to the constituent mass M in the homogeneous phase. 
Moreover, in the homogeneous chirally broken phase, M remains equal to the vacuum mass M^ac as 
long as fi < M^ac- Hence, if at T = there is a second-order phase transition to the inhomogeneous 
phase at fic,inh < Myac, then a single-soliton solution exists at /i = fic,inh, with the density profile 
given by Eq. ([19]) and A = Myac- 

In the next section this will be investigated further from the energy-per-particle perspective. 



III. ENERGY PER PARTICLE 

Starting from the thermodynamic potential, other thermodynamic quantities can be derived 
in the usual way, as long as we are only interested in spatial averages. Restricting ourselves to 
zero temperature, this means that the pressure p, the averaged quark number density n and the 
averaged energy density e are given by 

dp 

p = -Qyac) , e = -p + pn. (20) 

Here VL^ac is the value of the thermodynamic potential at its minimum in vacuum, which we 
subtract to define the vacuum pressure to be zero. As a consequence, the energy density of the 
vacuum vanishes as well. The average energy per quark is then given by 

1 = 4 = -? + .. (21) 

An n 

In the context of the interpretation of quark droplets as "baryons" , the thermodynamics is often 
discussed in terms of the baryon number density pB = n/Nc and the energy per baryon E/A = 
NcE/N, see, e.g., Refs. BSH- 

However, for most quantities we are going to discuss in this 
article, it is more natural to work with quark number densities and E/N. We therefore keep the 
notation introduced above, noting that the conversion to baryon quantities is simply a factor of 
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FIG. 2. Pressure as function of the chemical potential (left) and energy per quark as a function of the 
spatially averaged quark number density (right). The homogeneous chirally broken and restored solutions 
are indicated by the solid and dotted lines, respectively, while the dashed lines indicate the inhomogeneous 
solitonic solutions. 

Nc = 3. Moreover, in our numerical examples we will scale the densities by uq = NcPo, so that 
n/no = Pb/po- Here po = 0.17 fm~^ is the nuclear matter saturation density, i.e., the corresponding 
quark number density is no = 0.51 fm~'^. 

From Eq. (|2ip . it follows that the density derivative of E/N is given by 

where we have used that p = de/dn at T = and fixed volume. For n 7^ 0, this means that E/N 
has an extremum at the points where the pressure vanishes, and it takes the value E/N = p at 
these points. 

For n — )• 0, on the other hand, the exact behavior of E/N depends on the density dependence 
of the pressure. In the case of homogeneous quark matter, the NJL model at low densities behaves 
like an ideal nonrelativistic gas of constituent quarks, p oc n^/^. Consequently, E/N goes to p, 
which in turn converges to the vacuum constituent quark mass Myac, while the density derivative 
of E/N diverges at n = 0. As we will see below, the behavior of inhomogeneous matter is rather 
different. 

To this end, we now turn to the numerical results. Our model contains two parameters: the 
coupling constant G and the Pauli-Villars regulator A. We fix them by fitting the pion decay 
constant in vacuum to its value in the chiral limit, f^^ = 88 MeV, and by choosing a reasonable value 
for the constituent quark mass in vacuum. If not stated otherwise, we choose M^ac = 330 MeV, 
corresponding to A = 728.368 MeV and GA^ = 6.599. 

In the left panel of Fig. [21 we show the pressure as a function of p. The homogeneous chirally 
broken solutions are indicated by the solid line, where the upper branch corresponds to the minima 
of the thermodynamic potential, i.e., to the stable or metastable solutions, while the lower branch 
corresponds to the maxima, i.e., to the unstable solutions. The chirally restored solutions are 
indicated by the dotted line. Restricting the analysis to these homogeneous solutions, we find 
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a first-order chiral phase transition at /i — l^c,hom — 

329.9 MeV, i.e., slightly below n = Myac- 
Accordingly, the energy per particle in the restored phase, indicated by the dotted line in the 
right panel of Fig. [21 has a minimum with E/N = ^c,hom at n = (2/7r^)/x^ whereas in the 
homogeneous chirally broken solution (solid line), E/N is always larger, converging to Myac at 
n = with an infinite slope. Thus, in the "old picture", we would interpret the minimum in the 
restored phase as a bag-model-like quark droplet with a binding energy per quark of M^^q — tJ'c,hom' 
This picture is changed if we allow for the one-dimensional solitonic solutions, as indicated by 
the dashed lines in Fig. [2j We then find a second-order phase transitiorlfl from the homogeneous 
chirally broken phase to the inhomogeneous phase at /i = ^c,inh = 323.2 MeV (left panel). As 
discussed in Sect. Ill (J\ the inhomogeneous phase at this point corresponds to a single soliton 
{v = 1) with vanishing background density. Hence, there is no longer a stable solution with zero 
pressure and nonzero average density, and therefore the only minimum of the energy per particle 
exists at n = (right panel). On the other hand, the localized quarks inside the soliton experience 
additional binding, so that E/N does not go to M^ac at n = 0, as for homogeneous matter, but to 
Mc,m/i) which is smaller than Myac ™ this example. The binding effect is also visible at nonvanishing 
n. In particular, the chirally restored solution with the minimal E/N is unstable against forming 
a soliton lattice with the same average density. We find that here the solitons still have a sizable 
overlap, with density peaks separated by about a = 1.5 fm. This system can then lower its energy 
further by expansion. 

In this context, a striking difference to the homogeneous case is the fact that the density 
derivative of E/N does not diverge at n = but, on the contrary, the function is extremely flat. 
According to Eq. (I22p . this means that the pressure goes to zero with a high power of n. Further 
insight can be obtained from the observation in Ref. [IS] that the density rise above the onset of 
the solitonic phase is consistent with the parametrization 

n(/^) = -rT^^^^^^^' (23) 
ln(///^c,mfc - 1) 

where c is a constant parameter. This formula was motivated by a similar behavior in the Gross- 
Neveu model. Strictly speaking, it describes the density change 5n relative to the density at 
A* = fJ'c,inh- However, since Hcinh < Myac in the present case, we have n[^c,inh) = and, hence, 
5fi = n. We then find 

dp dp f dn\~'^ C/i^,m/i -c/x3 ,/n /r,.N 

IT = IT [it I = — = — ^ ' ' (24) 
on Ofi \Ofi J n 

where we have used Eqs. ([2D|) and ([25]) to evaluate the derivatives. It follows that ^ is exponentially 
suppressed for n — )• 0. The same is true for all higher derivatives and all derivatives of E/N, thus 
explaining its flatness. In fact, integrating Eq. (j24p to obtain |?(n) and inverting Eq. (j23p for fi{n), 



^ This means that the elhptic modulus v decreases continuously from v — 1 at fj,c,inh to smaller values inside the 
inhomogeneous phase. Obviously, this is impossible to prove by numerical calculations. Strictly speaking, we find 
that u does not drop discontinuously from 1 to a value smaller than 1 — 10~^^. At this point, the distance between 
the solitons is about a — 11 fm, which is well above the size of the solitons (see Fig. 



10 




FIG. 3. Energy per particle as a function of the average density at T = 0: numerical results (red solid line) 
and according to Eq. ([SS]) with c — 0.97985 and Hc.mh = 323.223 MeV (blue dashed line). 

we get from Eq. ([2T]) that the energy per particle at low densities should be given by 

/ . CO \ 

E _ 

f^c,inh 



3 

-I , -c/i3 . Jn ^^c,m/t /■ , e~ 
1 + e '^c.infi/ _ / (^x — 



n J X 
3 ' 



(25) 



In Fig. [3l this expression is compared with the numerical results for E/N. Fitting the parameters 
c and Hc,inh to the data below n = uq (left), we find a reasonable description up n = 2no (right), 
where the increase of E/N is more than a factor of 50 larger. We remark that the fitted value for 
c is very close to fJ-cinh/^vac, but we have not been able to show this analytically. 

For completeness, we also comment on the behavior at high densities. In our example with 
^vac = 330 MeV, the system stays inhomogeneous up to arbitrarily high chemical potentials. For 
somewhat lower values of M^ac, there is first a second-order phase transition from the solitonic 
phase to the restored phase, but the system gets inhomogeneous again at higher chemical poten- 



tials. As discussed in detail in Ref. 2l|, this so-called "inhomogeneous continent" is not a trivial 
regularization effect, but it cannot be excluded that it is a model artifact. In this article, however, 
we are mainly interested in the low-density behavior of the model, where this issue is irrelevant. 

The results shown in Fig. [2] have been obtained for specific parameters, and one might wonder 
how robust they are when these are changed. As discussed in Refs. 0, QjS]; the binding energy of 
homogeneous quark matter depends on the amount of scalar attraction, which can be parametrized 
by M^ac- In Fig- HI we therefore show the binding energies per quark in homogeneous matter and 
in the solitons, Ei,^hom = M^ac - f^cMm. and Eb^soi = M^ac - fJ-c,inh, respectively, as functions of 

Myac- 

We see that both curves start at the same point around M^ac = 250 MeV with a negative 
binding energy. For homogeneous matter, this point corresponds to the limiting case where the 
phase transition turns from first to second order when M^ac is lowered further. In other words, this 
point corresponds to the case where the tricritical point of the phase boundary in the T — fi plane 
is just located at the T = axis. Since in this model the tricritical point is equal to the Lifshitz 
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FIG. 4. Binding energy per quark for homogeneous matter (solid line) and for solitons (dashed line) as 
functions of the vacuum constituent quark mass. 



point [22|, i.e., the point where the two homogeneous phases and the inhomogeneous phase meet, 
the binding energies of homogeneous and inhomogeneous matter are equal at this point, and both 
solutions cease to exist at lower values of M„ac- 

When Myac is increased, the binding energies rise. The would-be first-order phase transition 
from the homogeneous chirally broken to the restored phase is now inside the inhomogeneous 
regime, i.e., fJ-c,inh < tJ'c,hom and, hence, Eh^soi > -E-s./iom- This means, the chirally restored solution 
with the lowest E/N is always unstable against forming a soliton lattice. 

On the other hand, for Myac < 315 MeV, 

f^c,inh is still smaller than Myac- Then the density 
in the homogeneous chirally broken phase is already nonzero when the phase transition to the 
inhomogeneous phase takes place. Thus, as discussed in Sec. Ill C[ the soliton is embedded in a 
homogeneous background of constituent quarks at this point. As the pressure is nonzero, the system 
wants to expand. However, since all solutions with a lower average density are homogeneous, this 
means that the inhomogeneous phase, including the single-soliton solution is not stable without 
applying external forcesQ As a consequence, the lowest E/N is obtained for a dilute gas of 
constituent quarks in the limit of zero density. 

For Myac > 315 MeV, Eh^soi is positive, i.e., the lowest E/N corresponds to a single soliton 
state, as discussed above. For Myac > 330 MeV, E^j^hom is positive as well. So there would be 
stable droplets of homogeneous matter if we could ignore the possibility of inhomogeneity. As 
explained above, this is, however, not the case. At Myac = 400 MeV, for example, the binding 
energy for homogeneous matter is about 30 MeV per quark, while it is about 40 MeV per quark 
for the solitons. 

For simplicity, all calculations in this article are done in the chiral limit. The mass functions 
and phase diagrams for non- vanishing bare quark masses have been investigated in Refs. [isl. [l6|. 

■* We have seen that inhomogeneous solutions with v = 1 and A = M are thermodynamically degenerate with 
homogeneous matter with mass M. Hence, one may argue that the single solitons survive also below ^ = ^c,inh, 
down to jj. — Myac, where they would have vanishing binding energy. On the other hand, they are no longer solutions 
of the gap equation §^ = 0, and it is therefore unclear whether they are self-consistent and thermodynamically 
consistent solutions. Here we choose not to further investigate this question, since in this work we are mainly 
interested in solutions with a positive binding energy, ^lc,inh < M^ac- 
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and turned out not to be very different. In particular, the inhomogeneous phase is deUmited by 
second-order phase boundaries, and the mass function at the boundary towards lower /j, takes the 
form of a single soliton. Therefore, we do not expect our results to change qualitatively if finite 
bare quark masses are considered. 



IV. PROPERTIES OF SINGLE SOLITONS 



After having explored the conditions for the existence of single self-bound solitons, we would 
now like to investigate their properties in more details. 

As discussed in Sec. Ill C\ the density profile is given by Eq. (fT9]l with fi = fic,inh and A = M^ac, 

= A^/A^c M,ac^^l,nh ^26) 
Svr cosh^ (M^iac 

and i^cinh < Myac- We have already noted that nsoi{z) decreases exponentially at large \z\ and 
therefore the average density n vanishes. On the other hand, the central density at z = is larger 
than the density of restored quark matter at the same chemical potential, 

nsol{0) = M^acfJ-cinh = ^rest , (27) 

where rirest = /^c ink- This can be interpreted as a bag-pressure effect, which pushes the 

quarks out of the chirally broken vacuum and squeezes them into the restored regions [16l | . 

The number of quarks in the soliton per transverse area A_i_ is obtained by integrating Eq. (j26|) 
over z. One finds 

17 - ^c,nh ■ (28) 

Another interesting quantity is the longitudinal rms "radius" , 

which is a measure for the half-size of the soliton in z-direction. 

Similarly, we can define the "soliton averaged density", i.e., the density- weighted integral over 
the density divided by the number of quarks, 

_ / dzn^,^i{z) _ NfNc 2 

= Jdznsoliz) = -T^^^-'^^'c,^nh■ (30) 

Hence 

{ri)sol = ^UsoliO) = j^^^Urest, (31) 

i.e., while the maximal density in a self-bound soliton is always larger than the density in the 
restored phase at the same chemical potential, {n)soi can be smaller. For instance, for Myac = 
400 MeV, we have Myac/l^c,inh = 1-12 < A/tt. 
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In order to interpret these expressions, we evaluate them exphcitly and compare the results with 
the droplet picture for homogeneous matter. As discussed earlier, the most stable homogeneous 
solution corresponds to quark matter in the restored phase at the critical chemical potential iJ-c,hom, 
provided fic,hom < ^vac- The density is, thus, given by 

_ ^f^<^ 3 

Assuming that the homogeneous solutions could be taken over to describe small quark matter 
droplets, the volume of a "baryon" with Nc quarks would be Vhom = ^c/nhom- For spherical bags, 
this would correspond to a radius of 

97r \ 1 



- ^ . (33) 

Similarly, we may restrict the transverse size of the soliton in such a way, that it contains 
quarks. Taking a circular shape, Eq. (p8|) yields 



Rsoi,± = ^= ^— . (34) 

Of course, this is a different geometry than the spherical shape we have assumed in the homogeneous 
case. To achieve a better comparability, it is therefore better to consider cylindrical droplets of 
homogeneous matter. For simplicity, we assume that the cylinder has a transverse radius Rhom,c 
and a height 2Rhom,c- Since the volume must remain unchanged, Rhom,c is then related to the 
radius of the sphere by Rhom,c = (2/3)^/^i?fcom,s, i-e., 

3vr \ 1 



Rhom,c - lyrr ■ (35) 

Moreover, for a better comparability with the longitudinal rms radius -R^J^V translate 
Rsoi,± and Rhom,c into rms radii. This yields 

'Rsol,± ~ Rsol,± ) Rhom,± ~ "^J^ Rhom,c i -^ftom, || ~ ~^J% Rhom,c : (36) 

where the different prefactors for longitudinal and transverse radii are due to the different dimen- 
sionality of these directions. 

Finally, we define a "baryon mass" for the homogeneous and solitonic solutions as Nc times 
E/N, i.e., 

Mb, horn = NcHchom > Mb, sol = NcHcinh ■ (37) 

The results for Nj = 2 and A'^c = 3 and two different values of M^ac are summarized in Table [J 
for homogeneous droplets and in Table HH for the solitons. We find that the qualitative and even the 
quantitative behavior is similar for both cases. The "baryon masses" rise with the vacuum quark 
mass, but are below 3Myac because of binding effects. Since the binding increases with increasing 
quark masses, the densities increase as well, while the radii decrease. As discussed in the previous 
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Myac [MeV] 


Mc./iom [MeV] 


MBJio,n [MeV] Uhom/no Rhom,s [fm] Rhol,± 


[fm] i?— [fm] 




330 

400 


329.9 
371.0 


yoy.i x.ou u.yi u.^ju 

1113.0 2.64 0.81 0.50 


0.41 




TABLE I. Properties of homoj 
mass Myac- 


;eneous "baryon droplets" for two values of the vacuum 


constituent quark 




Myac [MeV] 


Aic.in/i [MeV] 


Mb^soI [MeV] nsoi{0)/na {n)soi/no Rsoi,± [fm] 


Rsol,± i™[ RsoLW 


[fm] 


330 
400 


323.2 
357.4 


969.7 2.10 1.40 0.86 
1072.2 3.11 2.08 0.78 


0.61 0.54 
0.55 0.45 



TABLE II. Properties of solitonic "baryons" for two values of the vacuum constituent quark mass M, 



section, the solitons are bound more strongly than homogeneous matter. As a consequence, the 
solitons have a larger central density, despite the fact that the chemical potential is lower. The 
soliton averaged density {n)soU on the other hand, is smaller than the density in homogeneous 
droplets. Therefore, the solitons have larger rms radii. 

Nevertheless, the general agreement of the various rms radii for a given quark mass turns out 
to be quite good. In particular, it is remarkable that , which is an intrinsic property of the 
soliton, comes out to be similar to the other radii, which have been introduced by hand in order 
to have three quarks in a "baryon". On the other hand, the sharp radii Rhom,s and Rsoi,± are 
considerably larger, showing that the numbers are rather sensitive to the used definition of the 
radius. 



V. SOLITON-ANTISOLITON INTERACTIONS 

Having discussed the properties of single solitons, we now move away from this limit and inves- 
tigate what happens when the solitons approach each other. 

As explained in Sec. [Ill the inhomogeneous solutions are characterized by the parameters u and 
A, which are obtained by minimizing the thermodynamic potential at given T and fi. In particular 
the distance o between the neighboring solitons and antisolitons depends on u and A, as detailed 
in Eq. (jl2p . This allows us to plot the thermodynamic quantities of the system as functions of a, 
which is sometimes more instructive than plotting them against or n. 

As before, we limit ourselves to T = 0. At the boundary to the homogeneous chirally broken 
phase, we have v = \, corresponding to a — )• oo, while with increasing chemical potential the 
distance quickly becomes smaller. For large distances, the density distribution of the soliton lattice 
does not differ much from a linear superposition of single solitons. The average density is therefore 
given by 

2 

noo{a) = —. , (38) 

47r a 

i.e., the column density of a single soliton, Eq. (j28p . divided by the distance. At smaller distances, 
on the other hand, the interaction between the solitons leads to nonlinearities, giving rise to 
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FIG. 5. Left: Ratio of the average density n and the the corresponding value fioo, Eq- ([38|). for noninteracting 
sohtons as a function of the sohton-antisohton distance a. Right: Energy per nucleon (sohd hue) and force 
per nucleon (dashed line) as functions of a. 



deviations from the trivial 1/a behavior. This is shown in the left panel of Fig. [5l where n/fioo is 
displayed as a function of a. We see that the ratio is very close to unity for a > 2 fm and rises 
steeply when the distance is decreased below 1 fm. Here a is smaller than 2i?^^||, i.e., the solitons 
strongly overlap. 

In order to describe the soliton-antisoliton interactions more quantitatively, we convert the 
density axis of Fig. [2] (left) into the corresponding distance and plot E/N as a function of a. The 
result is displayed in the right panel of Fig. [5j We also show its derivative with respect to the 
distance, which corresponds the force per quark by which the solitons repel each other. 

According to the flatness of E/N at low average densities, we now find a flat behavior and, hence, 
an almost vanishing force at large distances, where "large" means again more than about 2 fm. At 
lower distances, the repulsion between the solitons becomes increasingly important. At a = 1 fm, 
for instance, we find a value of about 100 MeV/fm per quark. 



VI. INCLUDING VECTOR INTERACTIONS 



It is also interesting to study the influence of vector interactions, which are very important at 
finite density, as known, e.g., from the Walecka model 23|]. In the NJL model with homogeneous 
condensates, vector interactions have been shown to weaken the first-order chiral phase transition, 
and already at rather small values of the vector coupling, the phase transition turns into second 



order or a crossover 



24 



261 ]. In terms of E/N, this is easily understood from the fact that the 



vector interaction, described by a term 



(40) 
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in the Lagrangian, adds a term Gyn^ to the energy density, i.e., E/N is enhanced by Gyn [Tl]. 
Hence, the minimum in the restored phase at finite density gets increasingly disfavored with in- 
creasing Gv, whereas the energy at n = stays unaff'ected (see Ref. '2()] for a recent general 
discussion of this point). 

The effect of vector interactions on inhomogeneous phases has been investigated in Ref. p/7|]- In 
that analysis the approximation was made to replace the density n{z) in the mean-field Lagrangian 
by the spatial average n. This is a good approximation close to the restored phase and in particular 
at the Lifshitz point. The shape of the mass function at a given density is then independent of Gy, 
and the known analytical solutions for Gy = could basically be taken over. If we could apply 
the same approximation to our present analysis, we would obtain 



E 
N 



E 

N 



Gv=0 



+ Gvn^ {§] , (41) 



similar to the homogeneous case. This would further stabilize the minimum at n = 0. 

It is obvious, however, that the replacement of n by n is not a good approximation at low average 
densities where the quarks are strongly localized in the solitons and therefore feel a much stronger 
repulsion than suggested by Gyn- For instance, the energy of a single soliton is still enhanced 
by the vector repulsion, even when the homogeneous background density and, thus, the average 
density of the system goes to zero. Thus, the correction to E/N should rather be obtained by 
integrating the local correction to the energy density, 6e{z) = Gyn^(z) over the volume and divide 
it by the integrated quark number density. Since the integrals over the transverse area cancel, one 
obtains 



E 

N 



E 

N 



^ fdzn^iz) ( E\ 



This is still an approximation, at least as long the density profiles n{z) for Gy = which were 
given in Sec. Ill Bl are used. We expect that, in a fully self-consistent treatment, the vector repulsion 
between the quarks leads to a broadening of the density distribution, which lowers the energy. 
Eq. (j42p with the unmodified density profiles should therefore be taken as an upper limit of E/N, 
while Eq. (|4ip provides a lower limit. 

Making use of the periodicity of the soliton lattice, Eq. (j42p can be simplified to 



E\ _ E 
NJ^ ~ N 



+ ^ / dzn\z). (43) 



Gv=o «^ 







where a = L/2 is the distance between the solitons, introduced in Eq. (|12p . For the single-soliton 
limit with fj.c,inh < ^vac, we have E/N\gy=o = fJ-cinh, while the integrals in Eq. ()i2]) have readily 
been worked out in Sec. IIVI This yields 



(44) 



with {n)soi given in Eq. ([30 
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FIG. 6. Energy per particle as a function of the average density for a vector coupling Gy — G/2. The 
shaded area marks the range between the upper limit Eq. (j43|) and the lower limit Eq. (|4T|). The n — 0- value 
of the upper limit, Eq. (j44]), is denoted by the dash-dotted line. Also shown are E/N for homogeneous 
matter (solid line) and for inhomogeneous matter at Gy = (dashed line). 



In Fig. [6] our results for Gy = G/2 are displayed as functions of the average density. The range 
between the upper and lower limits of E/N is indicated by the shaded area. For comparison we 
also show the results for homogeneous matter and for inhomogeneous matter at Gy = 0. One can 
see that at high densities (E/N)^, {E/N)<^, and E/N for homogeneous matter become practically 
degenerate. This is not surprising, since in this regime the amplitude A of the mass function 
becomes small and the density profile gets more and more washed out 16(]. At intermediate densities 
we find the energy of homogeneous matter to be higher than the upper limit of inhomogeneous 
matter, i.e., the inhomogeneous solution should be favored in this region. 



The situation is less clear at lower densities. In the zero-density limit, (E/N)^ and [E/N) 



hom 



converge against the corresponding limits without vector interactions, i.e., fic ink 

and Myac, respec- 



tively, while (E/N)^ approaches the value given in Eq. 
small, 



. If the vector coupling is sufficiently 



Gv < 



vac fJ'cAnh 



(n) 



sol 



NfNcflc,inh \fJ'C,inh M, 



(45) 



the energy of homogeneous matter remains above the upper limit for inhomogeneous matter, even 
at n = 0, i.e., we can be rather sure that the inhomogeneous solutions stay favored. In the present 
example, however, Gy is not so small, and we find the ordering (E/N)^ < [E/N)hom < [E/N)y 
at low densities. If the correct inhomogeneous solution is close to the upper limit, this could mean 
that the ground state at low densities is homogeneous. On the other hand, it is also possible that 
the inhomogeneous solution remains favored if the solitons change their size in reaction to the 
repulsive vector interaction. 
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VII. CONCLUSIONS 



In this article, we have studied the existence and the properties of self-bound quark matter in 
the NJL model at zero temperature, focusing on inhomogeneous structures with one-dimensional 
spatial modulations. The analysis was done in mean-field approximation. 

For homogeneous matter, it was found long time ago that the model seems to allow for stable 
"droplets" of quark matter in the chirally restored phase if the interaction is sufficiently attractive. 
These droplets have vanishing pressure and a chemical potential lower than the vacuum constituent 
quark mass, so that they are in mechanical and chemical equilibrium with the vacuum. Related to 
this, they correspond to a minimum of the energy per particle as a function of density, so that they 
are stable against homogeneous expansion or collapse. Neglecting finite size effects, this suggests 
to interpret these solutions as quark bags, and the natural expectation would be that they have a 
spherical shape if surface effects are taken into account. 

Allowing for one-dimensional inhomogeneities, however, it turns out that the homogeneous 
droplets are unstable against forming a soliton lattice. The solitons repel each other, so that 
the state with the lowest energy per particle is reached at infinite lattice spacing, corresponding 
to a vanishing spatially averaged density. Inside the solitons, on the other hand, the density is 
finite, roughly of the same order as in the homogeneous droplets. Their longitudinal size is about 
1 fm, determined by the inverse of the vacuum constituent quark mass. Being one-dimensional 
objects embedded in the three-dimensional space, the solitons are infinite in the two transverse 
directions. Thus, taking these results as they are, quark matter at low average density should have 
a lasagne-like structure, with parallel plates of high densities and voids in between. 

At this point, we should ask ourselves how these results can be interpreted. In QCD, we 
expect that compressed quark matter, when it is released, will expand and finally hadronize. At 
zero temperature, this means that the matter should split up into baryons, each consisting of Nc 
valence quarks. These baryons may further interact with each other, forming nuclei or nuclear 
matter, but keep their individuality as separate color-singlet objects. 

In the NJL model, the "droplet" solutions found in the analysis of homogeneous quark matter 
have been suggested to be interpreted as schematic baryons, since they are stabilized by the bag 
pressure and have a reasonable density. Of course, strictly speaking, these solutions are infinite 
objects, and a separation into finite baryons would require a negative surface tension Q], while 
recent analyses suggest that it is positive 1CH12||. Prom this perspective, the one-dimensional 
solitons look like a step in the right direction, as they are at least finite in one dimension, where 
they have a reasonable size. In particular, one might hope that the consideration of higher- 
dimensional inhomogeneous phases could reveal further instabilities, eventually leading to finite 
localized baryons as the true ground state of matter at low densities. 

Unfortunately, this does not seem to be the case. Phases with two-dimensional modulations 
have been studied in Ref. 17| and were found to be disfavored against one-dimensional modulations 



at low densities. Although the analysis was restricted to sinusoidal shapes and certain parameters, 
it is unlikely that this will change if other shapes or parameters, or even three-dimensional modu- 
lations are considered. Nevertheless, more systematic studies in this direction are highly desirable, 



19 



in particular since at nonzero ternperature one-dimensional periodic structures are known to be 



soliton model 



unstable against fluctuations 27|, |28[ . One should also revisit the old works on the chiral quark 



29l43ll| and work out their relation to the present model. 



Of course, there is a priori no reason to expect finite baryons to be the most favored objects 
in a nonconfining model. On the other hand, the model predictions may still have some relevance 
in the deconfined phase. The emergence of one-dimensional modulations can be understood as a 



relic of the Peierls instability in 1 -|- 1 dimensions 84], which is a rather general mechanism. Also 
the fact that the longitudinal size and the internal density of the one-dimensional solitons are of 
the order to be expected for baryons might indicate that confinement effects are not very drastic. 
It is thus conceivable that lasagne-like patterns are preformed in expanding quark matter before 
hadronization takes place, and it would be interesting to work out possible observable signatures. 

The present calculations could also be improved in several aspects: In Sec. IVH we gave only a 
lower and an upper limit for the effect of vector interactions on the energy per particle. For the 
upper limit, which is probably closer to the true solution, we assumed that the density profiles 
remain unchanged when the vector interactions are switched on. However, we expect that the 
repulsive interaction leads to a broadening of the density peaks, which would lower the energy of 
the system. In this way the solitons may continuously go over into homogeneous matter, when 
the vector coupling is increased. We have also neglected the effect of spacelike vector condensates, 
which should be present in anisotropic systems. 

Moreover, we should allow for BCS pairing of the quarks in the solitons. Inside the solitons we 
find densities of two to three times nuclear-matter density, for which gaps of the order of 50 to 
100 MeV have been found in homogeneous quark matter. It would be interesting to see how this 
is changed for an inhomogeneous environment. 

Unfortunately, these improvements of the model can no longer be done by making use of the 
analytically known solutions of the 1-1-1 dimensional Gross-Neveu model, so that brute-force 
numerical diagonalizations of the Hamiltonian seem to be unavoidable. 

Finally, we should also include fluctuations. An interesting scenario would be that they leave 
the inhomogeneous phase (potentially with a higher-dimensional structure) intact but turn the 
second-order phase transition from the vacuum phase into first order. The minimum of E/N 
would then be shifted to nonvanishing average density. The resulting crystal could be a first step 
towards nuclear matter. 
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